Jeśli korzystasz z Google Colab skopiuj plik feature_names.json do katalogu głównego projektu.
Na tym laboratorium wykorzystamy zbiór danych Polish companies bankruptcy. Dotyczy on klasyfikacji, na podstawie danych z raportów finansowych, czy firma zbankrutuje w ciągu najbliższych kilku lat. Jest to zadanie szczególnie istotne dla banków, funduszy inwestycyjnych, firm ubezpieczeniowych itp., które z tego powodu zatrudniają licznie data scientistów. Zbiór zawiera 64 cechy, obliczone przez ekonomistów, którzy stworzyli ten zbiór, są one opisane na podlinkowanej wcześniej stronie. Dotyczą one zysków, posiadanych zasobów oraz długów firm.
Ściągnij i rozpakuj dane (Data Folder -> data.zip) do katalogu data obok tego notebooka. Znajduje się tam 5 plików w formacie .arff, wykorzystywanym głównie przez oprogramowanie Weka. Jest to program do "klikania" ML w interfejsie graficznym, jakiś czas temu popularny wśród mniej technicznych data scientistów. W Pythonie ładuje się je za pomocą bibliotek SciPy i Pandas.
Jeśli korzystasz z Linuksa możesz skorzystać z poniższych poleceń do pobrania i rozpakowania tych plików.
# !mkdir -p data
# !wget https://archive.ics.uci.edu/static/public/365/polish+companies+bankruptcy+data.zip -O data/data.zip
# !unzip data/data.zip -d data
W dalszej części laboratorium wykorzystamy plik 3year.arff, w którym na podstawie finansowych firmy po 3 latach monitorowania chcemy przewidywać, czy firma zbankrutuje w ciągu najbliższych 3 lat. Jest to dość realistyczny horyzont czasowy.
Dodatkowo w pliku feature_names.json znajdują się nazwy cech. Są bardzo długie, więc póki co nie będziemy z nich korzystać.
import json
import os
import numpy as np
import pandas as pd
from scipy.io import arff
data = arff.loadarff(os.path.join("data", "3year.arff"))
with open("feature_names.json") as file:
feature_names = np.array(json.load(file))
X = pd.DataFrame(data[0])
Przyjrzyjmy się teraz naszym danym.
X.head()
| Attr1 | Attr2 | Attr3 | Attr4 | Attr5 | Attr6 | Attr7 | Attr8 | Attr9 | Attr10 | ... | Attr56 | Attr57 | Attr58 | Attr59 | Attr60 | Attr61 | Attr62 | Attr63 | Attr64 | class | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0.174190 | 0.41299 | 0.14371 | 1.3480 | -28.9820 | 0.60383 | 0.219460 | 1.1225 | 1.1961 | 0.46359 | ... | 0.163960 | 0.375740 | 0.83604 | 0.000007 | 9.7145 | 6.2813 | 84.291 | 4.3303 | 4.0341 | b'0' |
| 1 | 0.146240 | 0.46038 | 0.28230 | 1.6294 | 2.5952 | 0.00000 | 0.171850 | 1.1721 | 1.6018 | 0.53962 | ... | 0.027516 | 0.271000 | 0.90108 | 0.000000 | 5.9882 | 4.1103 | 102.190 | 3.5716 | 5.9500 | b'0' |
| 2 | 0.000595 | 0.22612 | 0.48839 | 3.1599 | 84.8740 | 0.19114 | 0.004572 | 2.9881 | 1.0077 | 0.67566 | ... | 0.007639 | 0.000881 | 0.99236 | 0.000000 | 6.7742 | 3.7922 | 64.846 | 5.6287 | 4.4581 | b'0' |
| 3 | 0.024526 | 0.43236 | 0.27546 | 1.7833 | -10.1050 | 0.56944 | 0.024526 | 1.3057 | 1.0509 | 0.56453 | ... | 0.048398 | 0.043445 | 0.95160 | 0.142980 | 4.2286 | 5.0528 | 98.783 | 3.6950 | 3.4844 | b'0' |
| 4 | 0.188290 | 0.41504 | 0.34231 | 1.9279 | -58.2740 | 0.00000 | 0.233580 | 1.4094 | 1.3393 | 0.58496 | ... | 0.176480 | 0.321880 | 0.82635 | 0.073039 | 2.5912 | 7.0756 | 100.540 | 3.6303 | 4.6375 | b'0' |
5 rows × 65 columns
X.dtypes
Attr1 float64
Attr2 float64
Attr3 float64
Attr4 float64
Attr5 float64
...
Attr61 float64
Attr62 float64
Attr63 float64
Attr64 float64
class object
Length: 65, dtype: object
X.describe()
| Attr1 | Attr2 | Attr3 | Attr4 | Attr5 | Attr6 | Attr7 | Attr8 | Attr9 | Attr10 | ... | Attr55 | Attr56 | Attr57 | Attr58 | Attr59 | Attr60 | Attr61 | Attr62 | Attr63 | Attr64 | |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| count | 10503.000000 | 10503.000000 | 10503.000000 | 10485.000000 | 1.047800e+04 | 10503.000000 | 10503.000000 | 10489.000000 | 10500.000000 | 10503.000000 | ... | 1.050300e+04 | 10460.000000 | 10503.000000 | 10474.000000 | 10503.000000 | 9.911000e+03 | 10486.000000 | 1.046000e+04 | 10485.000000 | 10275.000000 |
| mean | 0.052844 | 0.619911 | 0.095490 | 9.980499 | -1.347662e+03 | -0.121159 | 0.065624 | 19.140113 | 1.819254 | 0.366093 | ... | 6.638549e+03 | -0.530082 | -0.014817 | 3.848794 | 1.429319 | 5.713363e+02 | 13.935361 | 1.355370e+02 | 9.095149 | 35.766800 |
| std | 0.647797 | 6.427041 | 6.420056 | 523.691951 | 1.185806e+05 | 6.970625 | 0.651152 | 717.756745 | 7.581659 | 6.428603 | ... | 5.989196e+04 | 55.978608 | 18.684047 | 190.201224 | 77.273270 | 3.715967e+04 | 83.704103 | 2.599116e+04 | 31.419096 | 428.298315 |
| min | -17.692000 | 0.000000 | -479.730000 | 0.002080 | -1.190300e+07 | -508.120000 | -17.692000 | -2.081800 | -1.215700 | -479.730000 | ... | -7.513800e+05 | -5691.700000 | -1667.300000 | -198.690000 | -172.070000 | 0.000000e+00 | -6.590300 | -2.336500e+06 | -0.000156 | -0.000102 |
| 25% | 0.000686 | 0.253955 | 0.017461 | 1.040100 | -5.207075e+01 | 0.000000 | 0.002118 | 0.431270 | 1.011275 | 0.297340 | ... | 1.462100e+01 | 0.005137 | 0.006796 | 0.875560 | 0.000000 | 5.533150e+00 | 4.486075 | 4.073700e+01 | 3.062800 | 2.023350 |
| 50% | 0.043034 | 0.464140 | 0.198560 | 1.605600 | 1.579300e+00 | 0.000000 | 0.050945 | 1.111000 | 1.199000 | 0.515500 | ... | 8.822900e+02 | 0.051765 | 0.106880 | 0.953060 | 0.002976 | 9.952100e+00 | 6.677300 | 7.066400e+01 | 5.139200 | 4.059300 |
| 75% | 0.123805 | 0.689330 | 0.419545 | 2.959500 | 5.608400e+01 | 0.072584 | 0.142275 | 2.857100 | 2.059100 | 0.725635 | ... | 4.348900e+03 | 0.130010 | 0.271310 | 0.995927 | 0.240320 | 2.093600e+01 | 10.587500 | 1.182200e+02 | 8.882600 | 9.682750 |
| max | 52.652000 | 480.730000 | 17.708000 | 53433.000000 | 6.854400e+05 | 45.533000 | 52.652000 | 53432.000000 | 740.440000 | 11.837000 | ... | 3.380500e+06 | 293.150000 | 552.640000 | 18118.000000 | 7617.300000 | 3.660200e+06 | 4470.400000 | 1.073500e+06 | 1974.500000 | 21499.000000 |
8 rows × 64 columns
feature_names
array(['net profit / total assets', 'total liabilities / total assets',
'working capital / total assets',
'current assets / short-term liabilities',
'[(cash + short-term securities + receivables - short-term liabilities) / (operating expenses - depreciation)] * 365',
'retained earnings / total assets', 'EBIT / total assets',
'book value of equity / total liabilities', 'sales / total assets',
'equity / total assets',
'(gross profit + extraordinary items + financial expenses) / total assets',
'gross profit / short-term liabilities',
'(gross profit + depreciation) / sales',
'(gross profit + interest) / total assets',
'(total liabilities * 365) / (gross profit + depreciation)',
'(gross profit + depreciation) / total liabilities',
'total assets / total liabilities', 'gross profit / total assets',
'gross profit / sales', '(inventory * 365) / sales',
'sales (n) / sales (n-1)',
'profit on operating activities / total assets',
'net profit / sales', 'gross profit (in 3 years) / total assets',
'(equity - share capital) / total assets',
'(net profit + depreciation) / total liabilities',
'profit on operating activities / financial expenses',
'working capital / fixed assets', 'logarithm of total assets',
'(total liabilities - cash) / sales',
'(gross profit + interest) / sales',
'(current liabilities * 365) / cost of products sold',
'operating expenses / short-term liabilities',
'operating expenses / total liabilities',
'profit on sales / total assets', 'total sales / total assets',
'constant capital / total assets', 'profit on sales / sales',
'(current assets - inventory - receivables) / short-term liabilities',
'total liabilities / ((profit on operating activities + depreciation) * (12/365))',
'profit on operating activities / sales',
'rotation receivables + inventory turnover in days',
'(receivables * 365) / sales', 'net profit / inventory',
'(current assets - inventory) / short-term liabilities',
'(inventory * 365) / cost of products sold',
'EBITDA (profit on operating activities - depreciation) / total assets',
'EBITDA (profit on operating activities - depreciation) / sales',
'current assets / total liabilities',
'short-term liabilities / total assets',
'(short-term liabilities * 365) / cost of products sold)',
'equity / fixed assets', 'constant capital / fixed assets',
'working capital', '(sales - cost of products sold) / sales',
'(current assets - inventory - short-term liabilities) / (sales - gross profit - depreciation)',
'total costs / total sales', 'long-term liabilities / equity',
'sales / inventory', 'sales / receivables',
'(short-term liabilities * 365) / sales',
'sales / short-term liabilities', 'sales / fixed assets'],
dtype='<U115')
DataFrame zawiera 64 atrybuty numeryczne o zróżnicowanych rozkładach wartości oraz kolumnę "class" typu bytes z klasami 0 i 1. Wiemy, że mamy do czynienia z klasyfikacją binarną - klasa 0 to brak bankructwa, klasa 1 to bankructwo w ciągu najbliższych 3 lat. Przyjrzyjmy się dokładniej naszym danym.
pd.Series, usuwając je z macierzy X. Przekonwertuj go na liczby całkowite.Uwaga: sugerowane jest użycie if w podpunkcie 1, żeby można było tę komórkę bezpiecznie odpalić kilka razy.
import matplotlib.pyplot as plt
# function for plotting class distribution
def plot_class_distribution(y):
counts = y.value_counts()
percentages = y.value_counts(normalize=True) * 100
ax = percentages.plot.bar()
ax.set_title("Classes Frequency")
ax.set_xlabel("Class")
ax.set_ylabel("Frequency [%]")
labels = [f"{count} ({perc:.2f}%)" for count, perc in zip(counts, percentages)]
ax.bar_label(ax.containers[0], labels=labels)
plt.show()
# your_code
if "class" in X.columns:
y = X.pop("class")
y = y.astype(int)
plot_class_distribution(y)
assert "class" not in X.columns
print("Solution is correct!")
Solution is correct!
Jak widać, klasa pozytywna jest w znacznej mniejszości, stanowi poniżej 5% zbioru. Taki problem nazywamy klasyfikacją niezbalansowaną (imbalanced classification). Mamy tu klasę dominującą (majority class) oraz klasę mniejszościową (minority class). Pechowo prawie zawsze interesuje nas ta druga, bo klasa większościowa jest trywialna. Przykładowo, 99% badanych jest zdrowych, a 1% ma niewykryty nowotwór - z oczywistych przyczyn chcemy wykrywać właśnie sytuację rzadką (problem diagnozy jako klasyfikacji jest zasadniczo zawsze niezbalansowany). W dalszej części laboratorium poznamy szereg konsekwencji tego zjawiska i metody na radzenie sobie z nim.
Mamy sporo cech, wszystkie numeryczne. Ciekawe, czy mają wartości brakujące, a jeśli tak, to ile. Policzymy to z pomocą biblioteki Pandas i metody .isna(). Domyślnie operuje ona na kolumnach, jak większość metod w Pandasie. Sumę wartości per kolumna zwróci nam metoda .sum(). Jeżeli podzielimy to przez liczbę wierszy len(X), to otrzymamy ułamek wartości brakujących w każdej kolumnie.
Pandas potrafi też stworzyć wykres, z pomocą funkcji np. .plot.hist() czy .plot.bar(). Przyjmują one opcje formatowania wykresu, z których korzysta pod spodem biblioteka matplotlib.
na_perc = X.isna().sum() / len(X)
na_perc.plot.bar(title="Fraction of missing values per column", figsize=(15, 5))
<Axes: title={'center': 'Fraction of missing values per column'}>
Jak widać, cecha 37 ma bardzo dużo wartości brakujących, podczas gdy pozostałe cechy mają raczej niewielką ich liczbę. W takiej sytuacji najlepiej usunąć tę cechę, a pozostałe wartości brakujące uzupełnić / imputować (impute). Typowo wykorzystuje się do tego wartość średnią lub medianę z danej kolumny. Ale uwaga - imputacji dokonuje się dopiero po podziale na zbiór treningowy i testowy! W przeciwnym wypadku wykorzystywalibyśmy dane ze zbioru testowego, co sztucznie zawyżyłoby wyniki. Jest to błąd metodologiczny - wyciek danych (data leakage).
Podział na zbiór treningowy i testowy to pierwszy moment, kiedy niezbalansowanie danych nam przeszkadza. Jeżeli zrobimy to czysto losowo, to są spore szanse, że w zbiorze testowym będzie tylko klasa negatywna - w końcu jest jej aż >95%. Dlatego wykorzystuje się próbkowanie ze stratyfikacją (stratified sampling), dzięki któremu proporcje klas w zbiorze przed podziałem oraz obu zbiorach po podziale są takie same.
"Attr37" ze zbioru danych.shuffle), ze stratyfikacją, wykorzystując funkcję train_test_split ze Scikit-learn'a.SimpleImputer.Uwaga:
if w podpunkcie 1random_state=0, aby wyniki były reprodukowalne (reproducible)stratify oczekuje wektora klas.fit()), a potem zastosować te nauczone wartości na obu podzbiorach (treningowym i testowym)# your_code
from sklearn.impute import SimpleImputer
from sklearn.model_selection import train_test_split
if "Attr37" in X:
X = X.drop("Attr37", axis="columns")
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=0.2, random_state=0, shuffle=True, stratify=y
)
mean_imputer = SimpleImputer(strategy="mean")
X_train = mean_imputer.fit_transform(X_train)
X_test = mean_imputer.transform(X_test)
import numpy as np
assert "Attr37" not in X.columns
assert not np.any(np.isnan(X_train))
assert not np.any(np.isnan(X_test))
print("Solution is correct!")
Solution is correct!
Zanim przejdzie się do modeli bardziej złożonych, trzeba najpierw wypróbować coś prostego, żeby mieć punkt odniesienia. Tworzy się dlatego modele bazowe (baselines).
W naszym przypadku będzie to drzewo decyzyjne (decision tree). Jest to drzewo binarne z decyzjami if-else, prowadzącymi do klasyfikacji danego przykładu w liściu. Każdy podział w drzewie to pytanie postaci "Czy wartość cechy X jest większa lub równa Y?". Trening takiego drzewa to prosty algorytm zachłanny, bardzo przypomina budowę zwykłego drzewa binarnego. W każdym węźle wykonujemy:
Taki algorytm wykonuje się rekurencyjnie, aż otrzymamy węzeł czysty (pure leaf), czyli taki, w którym są przykłady z tylko jednej klasy. Typowo wykorzystywaną funkcją jakości (kryterium podziału) jest entropia Shannona - im niższa entropia, tym bardziej jednolite są klasy w węźle (czyli wybieramy podział o najniższej entropii).
Powyższe wytłumaczenie algorytmu jest oczywiście nieformalne i dość skrótowe. Doskonałe tłumaczenie, z interaktywnymi wizualizacjami, dostępne jest tutaj. W formie filmów - tutaj oraz tutaj. Dla drzew do regresji - ten film.

Warto zauważyć, że taka konstrukcja prowadzi zawsze do overfittingu. Otrzymanie liści czystych oznacza, że mamy 100% dokładności na zbiorze treningowym, czyli perfekcyjnie przeuczony klasyfikator. W związku z tym nasze predykcje mają bardzo niski bias, ale bardzo dużą wariancję. Pomimo tego drzewa potrafią dać bardzo przyzwoite wyniki, a w celu ich poprawy można je regularyzować, aby mieć mniej "rozrośnięte" drzewo. Film dla zainteresowanych.
W tym wypadku AI to naprawdę tylko zbiór if'ów ;)
Mając wytrenowany klasyfikator, trzeba oczywiście sprawdzić, jak dobrze on sobie radzi. Tu natrafiamy na kolejny problem z klasyfikacją niezbalansowaną - zwykła celność (accuracy) na pewno nie zadziała! Typowo wykorzystuje się AUC, nazywane też AUROC (Area Under Receiver Operating Characteristic), bo metryka ta "widzi" i uwzględnia niezbalansowanie klas. Wymaga ona przekazania prawdopodobieństwa klasy pozytywnej, a nie tylko binarnej decyzji.
Bardzo dobre i bardziej szczegółowe wytłumaczenie, z interktywnymi wizualizacjami, można znaleć tutaj. Dla preferujących filmy - tutaj.
Co ważne, z definicji AUROC, trzeba tam użyć prawdopodobieństw klasy pozytywnej (klasy 1). W Scikit-learn'ie zwraca je metoda .predict_proba(), która w kolejnych kolumnach zwraca prawdopodobieństwa poszczególnych klas.
DecisionTreeClassifier). Użyj entropii jako kryterium podziału.roc_auc_score).Uwaga:
random_state=0plt.show() z Matplotlibafrom sklearn.metrics import roc_auc_score
# function for AUROC calculation
def get_auroc_score(model, X, y):
y_score = model.predict_proba(X)[:, 1]
auroc_score = roc_auc_score(y, y_score)
print(f"AUROC score: {auroc_score}")
return auroc_score
# your_code
from sklearn import tree
from sklearn.metrics import roc_auc_score
from sklearn.tree import DecisionTreeClassifier
model = DecisionTreeClassifier(random_state=0, criterion="entropy")
model.fit(X_train, y_train)
auroc = get_auroc_score(model, X_test, y_test)
AUROC score: 0.7266899766899767
fig = plt.figure(figsize=(200, 160))
_ = tree.plot_tree(model, filled=True, feature_names=X.columns, rounded=True)
assert auroc > 0.7
print("Solution is correct!")
Solution is correct!
// skomentuj tutaj
Moim zdaniem osiągnięty AUROC nie jest wysoki. Wytrenoway model ma "jakąś" moc klasyfikacyjną - AUROC > 0.5, więc nie jest to totalna losowość. Przejrzałem parę artykułów (1, 2, 3) i wyszło, że nie ma tak naprawdę czegoś takiego jak obiektywnie dobry AUROC, bo jest to zależne od sytuacji, ale często przyjmuje się, że ~0.7 to akceptowalny wynik - czyli na pewno da się lepiej.
Bardzo często wiele klasyfikatorów działających razem daje lepsze wyniki niż pojedynczy klasyfikator. Takie podejście nazywa się uczeniem zespołowym (ensemble learning). Istnieje wiele różnych podejść do tworzenia takich klasyfikatorów złożonych (ensemble classifiers).
Podstawową metodą jest bagging:
![]()
Typowo klasyfikatory bazowe są bardzo proste, żeby można było szybko wytrenować ich dużą liczbę. Prawie zawsze używa się do tego drzew decyzyjnych. Dla klasyfikacji uśrednienie wyników polega na głosowaniu - dla nowej próbki każdy klasyfikator bazowy ją klasyfikuje, sumuje się głosy na każdą klasę i zwraca najbardziej popularną decyzję.
Taki sposób ensemblingu zmniejsza wariancję klasyfikatora. Intuicyjnie, skoro coś uśredniamy, to siłą rzeczy będzie mniej rozrzucone, bo dużo ciężej będzie osiągnąć jakąś skrajność. Redukuje to też overfitting.
Lasy losowe (Random Forests) to ulepszenie baggingu. Zaobserwowano, że pomimo losowania próbek boostrapowych, w baggingu poszczególne drzewa są do siebie bardzo podobne (są skorelowane), używają podobnych cech ze zbioru. My natomiast chcemy zróżnicowania, żeby mieć niski bias - redukcją wariancji zajmuje się uśrednianie. Dlatego używa się metody losowej podprzestrzeni (random subspace method) - przy każdym podziale drzewa losuje się tylko pewien podzbiór cech, których możemy użyć do tego podziału. Typowo jest to pierwiastek kwadratowy z ogólnej liczby cech.
Zarówno bagging, jak i lasy losowe mają dodatkowo bardzo przyjemną własność - są mało czułe na hiperparametry, szczególnie na liczbę drzew. W praktyce wystarczy ustawić 500 czy 1000 drzew i będzie dobrze działać. Dalsze dostrajanie hiperparametrów może jeszcze trochę poprawić wyniki, ale nie tak bardzo, jak przy innych klasyfikatorach. Jest to zatem doskonały wybór domyślny, kiedy nie wiemy, jakiego klasyfikatora użyć.
Dodatkowo jest to problem embarassingly parallel - drzewa można trenować w 100% równolegle, dzięki czemu jest to dodatkowo wydajna obliczeniowo metoda.
Głębsze wytłumaczenie, z interaktywnymi wizualizacjami, można znaleźć tutaj. Dobrze tłumaczy je też ta seria filmów.
RandomForestClassifier). Użyj 500 drzew i entropii jako kryterium podziału.Uwaga: pamiętaj o ustawieniu random_state=0. Dla przyspieszenia ustaw n_jobs=-1 (użyje tylu procesów, ile masz dostępnych rdzeni procesora). Pamiętaj też o przekazaniu prawdopodobieństw do metryki AUROC.
# your_code
from sklearn.ensemble import RandomForestClassifier
model = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
model.fit(X_train, y_train)
auroc = get_auroc_score(model, X_test, y_test)
AUROC score: 0.8994111948657404
assert auroc > 0.85
print("Solution is correct!")
Solution is correct!
// skomentuj tutaj
Otrzymaliśmy zdecydowanie lepszy wynik niż przy drzewie decyzyjnym. Bazując na ogólnie przyjętych normach jest to dobry wynik.
Jak zobaczymy poniżej, wynik ten możemy jednak jeszcze ulepszyć!
W przypadku zbiorów niezbalansowanych można dokonać balansowania (balancing) zbioru. Są tutaj 2 metody:
Undersampling działa dobrze, kiedy niezbalansowanie jest niewielkie, a zbiór jest duży (możemy sobie pozwolić na usunięcie jego części). Oversampling typowo daje lepsze wyniki, istnieją dla niego bardzo efektywne algorytmy. W przypadku bardzo dużego niezbalansowania można zrobić oba.
Typowym algorytmem oversamplingu jest SMOTE (Synthetic Minority Oversampling TEchnique). Działa on następująco:
k najbliższych przykładów dla próbki, typowo k=5
Taka technika generuje przykłady bardzo podobne do prawdziwych, więc nie zaburza zbioru, a jednocześnie pomaga klasyfikatorom, bo "zagęszcza" przestrzeń, w której znajduje się klasa pozytywna.
Algorytm SMOTE, jego warianty i inne algorytmy dla problemów niezbalansowanych implementuje biblioteka Imbalanced-learn.
Użyj SMOTE do zbalansowania zbioru treningowego (nie używa się go na zbiorze testowym!) (klasa SMOTE). Wytrenuj drzewo decyzyjne oraz las losowy na zbalansowanym zbiorze, użyj tych samych argumentów co wcześniej. Pamiętaj o użyciu wszędzie stałego random_state=0 oraz przekazaniu prawdopodobieństw do AUROC. Skomentuj wynik.
Wartość ROC drzewa decyzyjnego przypisz do zmiennej tree_roc, a lasu do forest_roc.
# your_code
from imblearn.over_sampling import SMOTE
smote = SMOTE(random_state=0)
X_reb, y_reb = smote.fit_resample(X_train, y_train)
plot_class_distribution(y_reb)
print(f"Decison tree:")
decision_tree = DecisionTreeClassifier(random_state=0, criterion="entropy")
decision_tree.fit(X_reb, y_reb)
tree_roc = get_auroc_score(decision_tree, X_test, y_test)
print()
print(f"Random forest:")
random_forest = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest.fit(X_reb, y_reb)
forest_roc = get_auroc_score(random_forest, X_test, y_test)
Decison tree: AUROC score: 0.70995670995671 Random forest: AUROC score: 0.9047644274917003
assert 0.6 < tree_roc < 0.8
assert 0.8 < forest_roc < 0.95
// skomentuj tutaj
Wyniki są porównywalne do tych, które otrzymaliśmy poprzednio. AUROC można interpretować w ten sposób, że jest to prawdopodobieństwo, na to że spośród losowo wybranego TP i losowo wybranego TN, model zaklasyfikuje wyżej TP niż TN. Ta metryka nie hierarchizuje klas, przez co dobrze radzi sobie w przypadku ich niezbalansowania. Z tego powodu na zbiorze zrównoważonym, gdzie też nie uznaje jednej klasy za ważniejszej, radzi sobie podobnie.
W dalszej części laboratorium używaj zbioru po zastosowaniu SMOTE do treningu klasyfikatorów.
Lasy losowe są stosunkowo mało czułe na dobór hiperparametrów - i dobrze, bo mają ich dość dużo. Można zawsze jednak spróbować to zrobić, a w szczególności najważniejszy jest parametr max_features, oznaczający, ile cech losować przy każdym podziale drzewa. Typowo sprawdza się wartości z zakresu [0.1, 0.5].
W kwestii szybkości, kiedy dostrajamy hiperparametry, to mniej oczywiste jest, jakiego n_jobs użyć. Z jednej strony klasyfikator może być trenowany na wielu procesach, a z drugiej można trenować wiele klasyfikatorów na różnych zestawach hiperparametrów równolegle. Jeżeli nasz klasyfikator bardzo dobrze się uwspółbieżnia (jak Random Forest), to można dać mu nawet wszystkie rdzenie, a za to wypróbowywać kolejne zestawy hiperparametrów sekwencyjnie. Warto ustawić parametr verbose na 2 lub więcej, żeby dostać logi podczas długiego treningu i mierzyć czas wykonania. W praktyce ustawia się to metodą prób i błędów.
max_features:[0.1, 0.2, 0.3, 0.4, 0.5]scoring)max_features. Jest to atrybut wytrenowanego GridSearchCV.auroc.Uwaga:
random_state=0 i n_jobs# your_code
from sklearn.model_selection import GridSearchCV
parameters = {"max_features": [0.1, 0.2, 0.3, 0.4, 0.5]}
random_forest_estimator = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest_GS = GridSearchCV(
estimator=random_forest_estimator,
cv=5,
param_grid=parameters,
verbose=2,
scoring="roc_auc",
n_jobs=-1,
)
random_forest_GS.fit(X_train, y_train)
Fitting 5 folds for each of 5 candidates, totalling 25 fits
GridSearchCV(cv=5,
estimator=RandomForestClassifier(criterion='entropy',
n_estimators=500, n_jobs=-1,
random_state=0),
n_jobs=-1, param_grid={'max_features': [0.1, 0.2, 0.3, 0.4, 0.5]},
scoring='roc_auc', verbose=2)In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook. GridSearchCV(cv=5,
estimator=RandomForestClassifier(criterion='entropy',
n_estimators=500, n_jobs=-1,
random_state=0),
n_jobs=-1, param_grid={'max_features': [0.1, 0.2, 0.3, 0.4, 0.5]},
scoring='roc_auc', verbose=2)RandomForestClassifier(criterion='entropy', max_features=0.4, n_estimators=500,
n_jobs=-1, random_state=0)RandomForestClassifier(criterion='entropy', max_features=0.4, n_estimators=500,
n_jobs=-1, random_state=0)print(f'Optimal max_features value: {random_forest_GS.best_params_["max_features"]}')
auroc = get_auroc_score(random_forest_GS, X_test, y_test)
Optimal max_features value: 0.4 AUROC score: 0.9217020353383989
assert 0.9 <= auroc <= 0.95
print("Solution is correct!")
Solution is correct!
// skomentuj tutaj
Samo dostrajanie w moim przypadku nie zajęło zbyt dużo czasu - testowałem na dwóch maszynach i wyszło między 3-7 minut, więc akceptowalnie. Wynik się poprawił, więc jak najbardziej było warto.
W praktycznych zastosowaniach data scientist wedle własnego uznana, doświadczenia, dostępnego czasu i zasobów wybiera, czy dostrajać hiperparametry i w jak szerokim zakresie. Dla Random Forest na szczęście często może nie być znaczącej potrzeby, i za to go lubimy :)
Random Forest - podsumowanie
Drugą bardzo ważną grupą algorytmów ensemblingu jest boosting, też oparty o drzewa decyzyjne. O ile Random Forest trenował wszystkie klasyfikatory bazowe równolegle i je uśredniał, o tyle boosting robi to sekwencyjnie. Drzewa te uczą się na całym zbiorze, nie na próbkach boostrapowych. Idea jest następująca: trenujemy drzewo decyzyjne, radzi sobie przeciętnie i popełnia błędy na częsci przykładów treningowych. Dokładamy kolejne, ale znające błędy swojego poprzednika, dzięki czemu może to uwzględnić i je poprawić. W związku z tym "boostuje" się dzięki wiedzy od poprzednika. Dokładamy kolejne drzewa zgodnie z tą samą zasadą.
Jak uczyć się na błędach poprzednika? Jest to pewna funkcja kosztu (błędu), którą chcemy zminimalizować. Zakłada się jakąś jej konkretną postać, np. squared error dla regresji, albo logistic loss dla klasyfikacji. Później wykorzystuje się spadek wzdłuż gradientu (gradient descent), aby nauczyć się, w jakim kierunku powinny optymalizować kolejne drzewa, żeby zminimalizować błędy poprzednika. Jest to konkretnie gradient boosting, absolutnie najpopularniejsza forma boostingu, i jeden z najpopularniejszych i osiągających najlepsze wyniki algorytmów ML.
Tyle co do intuicji. Ogólny algorytm gradient boostingu jest trochę bardziej skomplikowany. Bardzo dobrze i krok po kroku tłumaczy go ta seria filmów na YT. Szczególnie ważne implementacje gradient boostingu to XGBoost (Extreme Gradient Boosting) oraz LightGBM (Light Gradient Boosting Machine). XGBoost był prawdziwym przełomem w ML, uzyskując doskonałe wyniki i bardzo dobrze się skalując - był wykorzystany w CERNie do wykrywania cząstki Higgsa w zbiorze z pomiarów LHC mającym 10 milionów próbek. Jego implementacja jest dość złożona, ale dobrze tłumaczy ją inna seria filmików na YT.

Obecnie najczęściej wykorzystuje się LightGBM. Został stworzony przez Microsoft na podstawie doświadczeń z XGBoostem. Został jeszcze bardziej ulepszony i przyspieszony, ale różnice są głównie implementacyjne. Różnice dobrze tłumaczy ta prezentacja z konferencji PyData oraz prezentacja Microsoftu. Dla zainteresowanych - praktyczne aspekty LightGBM.
LGBMClassifier). Przekaż importance_type="gain" - przyda nam się to za chwilę.Pamiętaj o random_state, n_jobs i prawdopodobieństwach dla AUROC.
# your_code
from lightgbm import LGBMClassifier
LGBM_clf = LGBMClassifier(importance_type="gain", random_state=0, n_jobs=-1)
LGBM_clf.fit(X_reb, y_reb)
auroc = get_auroc_score(LGBM_clf, X_test, y_test)
[LightGBM] [Info] Number of positive: 8006, number of negative: 8006 [LightGBM] [Info] Auto-choosing col-wise multi-threading, the overhead of testing was 0.003486 seconds. You can set `force_col_wise=true` to remove the overhead. [LightGBM] [Info] Total Bins 16065 [LightGBM] [Info] Number of data points in the train set: 16012, number of used features: 63 [LightGBM] [Info] [binary:BoostFromScore]: pavg=0.500000 -> initscore=0.000000 AUROC score: 0.9433748070111706
assert 0.9 <= auroc <= 0.97
print("Solution is correct!")
Solution is correct!
// skomentuj tutaj
LGBM osiągnał znakomity wynik, a ponadto, porównując z innymi klasyfikatorami, okazał się niesamowicie szybki, co tym bardziej działa na plus.
Boosting dzięki uczeniu na poprzednich drzewach redukuje nie tylko wariancję, ale też bias w błędzie, dzięki czemu może w wielu przypadkach osiągnąć lepsze rezultaty od lasu losowego. Do tego dzięki znakomitej implementacji LightGBM jest szybszy.
Boosting jest jednak o wiele bardziej czuły na hiperparametry niż Random Forest. W szczególności bardzo łatwo go przeuczyć, a większość hiperparametrów, których jest dużo, wiąże się z regularyzacją modelu. To, że teraz poszło nam lepiej z domyślnymi, jest rzadkim przypadkiem.
W związku z tym, że przestrzeń hiperparametrów jest duża, przeszukanie wszystkich kombinacji nie wchodzi w grę. Zamiast tego można wylosować zadaną liczbę zestawów hiperparametrów i tylko je sprawdzić - chociaż im więcej, tym lepsze wyniki powinniśmy dostać. Służy do tego RandomizedSearchCV. Co więcej, klasa ta potrafi próbkować rozkłady prawdopodobieństwa, a nie tylko sztywne listy wartości, co jest bardzo przydatne przy parametrach ciągłych.
Hiperparametry LightGBMa są dobrze opisane w oficjalnej dokumentacji: wersja krótsza i wersja dłuższa. Jest ich dużo, więc nie będziemy ich tutaj omawiać. Jeżeli chodzi o ich dostrajanie w praktyce, to przydatny jest oficjalny guide oraz dyskusje na Kaggle.
RandomizedSearchCV):param_grid = {
"n_estimators": [100, 250, 500],
"learning_rate": [0.05, 0.1, 0.2],
"num_leaves": [31, 48, 64],
"colsample_bytree": [0.8, 0.9, 1.0],
"subsample": [0.8, 0.9, 1.0],
}
classification_report), dla modelu LightGBM bez i z dostrajaniem hiperparametrów.auroc.Uwaga:
verbose=-1 przy tworzeniu LGBMClassifier, żeby uniknąć kolosalnej ilości logów, która potrafi też wyłączyć Jupyteraimportance_type, random_state=0 i n_jobs, oraz ewentualnie verbose w RandomizedSearchCV dla śledzenia przebiegun_jobs dla grid searcha będzie szybsze niż dla samego LightGBM; odpowiada to tuningowi wielu klasyfikatorów równolegle, przy wolniejszym treningu pojedynczych klasyfikatorówn_jobs=-1, bo wtedy stworzysz więcej procesów niż rdzeni i spowodujesz thread contentionfrom sklearn.metrics import classification_report
# function for printing the classification report
def get_classification_report(model, X_test, y_test, title):
y_pred = model.predict(X_test)
print(title)
print(classification_report(y_true=y_test, y_pred=y_pred))
get_classification_report(
LGBM_clf,
X_test,
y_test,
"Classification report for LGBM without hyper-parameter tuning:",
)
Classification report for LGBM without hyper-parameter tuning:
precision recall f1-score support
0 0.98 0.98 0.98 2002
1 0.60 0.60 0.60 99
accuracy 0.96 2101
macro avg 0.79 0.79 0.79 2101
weighted avg 0.96 0.96 0.96 2101
# your_code
from sklearn.model_selection import RandomizedSearchCV
param_grid = {
"n_estimators": [100, 250, 500],
"learning_rate": [0.05, 0.1, 0.2],
"num_leaves": [31, 48, 64],
"colsample_bytree": [0.8, 0.9, 1.0],
"subsample": [0.8, 0.9, 1.0],
}
LGBM_estimator = LGBMClassifier(
importance_type="gain", random_state=0, n_jobs=-1, verbose=-1
)
LGBM_RS = RandomizedSearchCV(
estimator=LGBM_estimator,
cv=5,
param_distributions=param_grid,
verbose=0,
scoring="roc_auc",
n_iter=30,
)
LGBM_RS.fit(X_reb, y_reb)
RandomizedSearchCV(cv=5,
estimator=LGBMClassifier(importance_type='gain', n_jobs=-1,
random_state=0, verbose=-1),
n_iter=30,
param_distributions={'colsample_bytree': [0.8, 0.9, 1.0],
'learning_rate': [0.05, 0.1, 0.2],
'n_estimators': [100, 250, 500],
'num_leaves': [31, 48, 64],
'subsample': [0.8, 0.9, 1.0]},
scoring='roc_auc')In a Jupyter environment, please rerun this cell to show the HTML representation or trust the notebook. RandomizedSearchCV(cv=5,
estimator=LGBMClassifier(importance_type='gain', n_jobs=-1,
random_state=0, verbose=-1),
n_iter=30,
param_distributions={'colsample_bytree': [0.8, 0.9, 1.0],
'learning_rate': [0.05, 0.1, 0.2],
'n_estimators': [100, 250, 500],
'num_leaves': [31, 48, 64],
'subsample': [0.8, 0.9, 1.0]},
scoring='roc_auc')LGBMClassifier(colsample_bytree=0.9, importance_type='gain', learning_rate=0.2,
n_estimators=500, n_jobs=-1, num_leaves=48, random_state=0,
verbose=-1)LGBMClassifier(colsample_bytree=0.9, importance_type='gain', learning_rate=0.2,
n_estimators=500, n_jobs=-1, num_leaves=48, random_state=0,
verbose=-1)print(f"Best parameters: {LGBM_RS.best_params_}\n")
get_classification_report(
LGBM_RS,
X_test,
y_test,
"Classification report for LGBM with hyper-parameter tuning:",
)
Best parameters: {'subsample': 1.0, 'num_leaves': 48, 'n_estimators': 500, 'learning_rate': 0.2, 'colsample_bytree': 0.9}
Classification report for LGBM with hyper-parameter tuning:
precision recall f1-score support
0 0.98 0.99 0.99 2002
1 0.80 0.54 0.64 99
accuracy 0.97 2101
macro avg 0.89 0.76 0.81 2101
weighted avg 0.97 0.97 0.97 2101
auroc = get_auroc_score(LGBM_RS, X_test, y_test)
assert 0.9 <= auroc <= 0.99
AUROC score: 0.9397925307016216
// skomentuj tutaj
To czy jest to pożądane zjawisko, to zależy od naszej taktyki. Załóżmy, że bardziej nas interesuje precyzja i chcemy być maksylmanie skuteczni w naszych przewidywaniach bankructwa, bo fałszywy alarm może prowadzić np. do zamknięcia dobrze prosperującej firmy albo drastycznej restrukturyzacji. W takim wypadku, dostrojenie hiper-parametrów miało duże znaczenie. Precyzja znacznie się polepszyła, co jest dla nas bardzo korzystne. Warto zwrócić uwagę na to, że AUROC w obu przypadkach jest dosyć podobny (~0.94). Podczas wykonywania tego labolatorium napotkałem masę różnych artykułów pt. "Why you should not use AUROC". Autorzy słusznie zwracali uwagę na potrzebę oceny modelu również innymi metrykami (często przewijał się Precision-Recall Curve), co w tym wypadku ma swoje odzwierciedlenie, jako że dla dwóch bliskich wyników AUROC, inne metryki, z perspektywy naszego biznesplanu - ważniejsze, mocno się różniły.
Boosting - podsumowanie
W ostatnich latach zaczęto zwracać coraz większą uwagę na wpływ sztucznej inteligencji na społeczeństwo, a na niektórych czołowych konferencjach ML nawet obowiązkowa jest sekcja "Social impact" w artykułach naukowych. Typowo im lepszy model, tym bardziej złożony, a najpopularniejsze modele boostingu są z natury skomplikowane. Kiedy mają podejmować krytyczne decyzje, to musimy wiedzieć, czemu predykcja jest taka, a nie inna. Jest to poddziedzina uczenia maszynowego - wyjaśnialna AI (explainable AI, XAI).
Taka informacja jest cenna, bo dzięki temu lepiej wiemy, co robi model. Jest to ważne z kilku powodów:
W szczególności można ją podzielić na globalną oraz lokalną interpretowalność (global / local interpretability). Ta pierwsza próbuje wyjaśnić, czemu ogólnie model działa tak, jak działa. Analizuje strukturę modelu oraz trendy w jego predykcjach, aby podsumować w prostszy sposób jego tok myślenia. Interpretowalność lokalna z kolei dotyczy predykcji dla konkretnych próbek - czemu dla danego przykładu model podejmuje dla niego taką, a nie inną decyzję o klasyfikacji.
W szczególności podstawowym sposobem interpretowalności jest ważność cech (feature importance). Wyznacza ona, jak ważne są poszczególne cechy:
Teraz będzie nas interesować globalna ważność cech. Dla modeli drzewiastych definiuje się ją bardzo prosto. Każdy podział w drzewie decyzyjnym wykorzystuje jakąś cechę, i redukuje z pomocą podziału funkcję kosztu (np. entropię) o określoną ilość. Dla drzewa decyzyjnego ważność to sumaryczna redukcja entropii, jaką udało się uzyskać za pomocą danej cechy. Dla lasów losowych i boostingu sumujemy te wartości dla wszystkich drzew. Alternatywnie można też użyć liczby splitów, w jakiej została użyta dana cecha, ale jest to mniej standardowe.
Warto zauważyć, że taka ważność cech jest względna:
Ze względu na powyższe, ważności cech normalizuje się często do zakresu [0, 1] dla łatwiejszego porównywania.
feature_names.Uwaga: Scikit-learn normalizuje ważności do zakresu [0, 1], natomiast LightGBM nie. Musisz to znormalizować samodzielnie, dzieląc przez sumę.
# function for plotting top 5 feature importances on a specific ax
def plot_feature_importances(feature_importances, ax, title):
head_indices = np.argsort(feature_importances)[-5:]
head = feature_importances[head_indices]
head_names = feature_names[head_indices]
ax.bar(head_names, head)
ax.set_ylabel("Feature Importance")
ax.set_title(title)
ax.tick_params(rotation=70, axis="x")
# your_code
fig, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(30, 9))
tree_feature_importances = decision_tree.feature_importances_
forest_feature_importances = random_forest.feature_importances_
LGBM_feature_importances = LGBM_clf.feature_importances_
LGBM_feature_importances /= LGBM_feature_importances.sum()
plot_feature_importances(
tree_feature_importances, ax1, "Feature Importance from Decision Tree [Top 5]"
)
plot_feature_importances(
forest_feature_importances, ax2, "Feature Importance from Random Forest [Top 5]"
)
plot_feature_importances(
LGBM_feature_importances, ax3, "Feature Importance from LGBM [Top 5]"
)
// skomentuj tutaj
W trzech różnych metodach zostały wyróżnione całkiem "logiczne" cechy: takie jak sales(n)/sales(n-1) oraz profit on operating activites/financial expenses. To ważne cechy, których odpowiednie stosunki (tzn. $\geq 1$) świadczą o dobrej kondycji firmy. Pierwsza kładzie nacisk na zdolność do rozwoju, a druga na wydajność finansową firmy (potrzebną np. w przypadku długów). Myslę, że wybór tych cech jak najbardziej ma sens.
Najpopularniejszym podejściem do interpretowalności lokalnych jest SHAP (SHapley Additive exPlanations), metoda oparta o kooperatywną teorię gier. Traktuje się cechy modelu jak zbiór graczy, podzielonych na dwie drużyny (koalicje): jedna chce zaklasyfikować próbkę jako negatywną, a druga jako pozytywną. O ostatecznej decyzji decyduje model, który wykorzystuje te wartości cech. Powstaje pytanie - w jakim stopniu wartości cech przyczyniły się do wyniku swojej drużyny? Można to obliczyć jako wartości Shapleya (Shapley values), które dla modeli ML oblicza algorytm SHAP. Ma on bardzo znaczące, udowodnione matematycznie zalety, a dodatkowo posiada wyjątkowo efektywną implementację dla modeli drzewiastych oraz dobre wizualizacje.
Bardzo intuicyjnie, na prostym przykładzie, SHAPa wyjaśnia pierwsza część tego artykułu. Dobrze i dość szczegółówo SHAPa wyjaśnia jego autor w tym filmie.
Wyjaśnialna AI - podsumowanie
Dokonaj selekcji cech, usuwając 20% najsłabszych cech. Może się tu przydać klasa SelectPercentile. Czy Random Forest i LightGBM (bez dostrajania hiperparametrów, dla uproszczenia) wytrenowane bez najsłabszych cech dają lepszy wynik (AUROC lub innej metryki)?
Wykorzystaj po 1 algorytmie z 3 grup algorytmów selekcji cech:
chi2 i mutual_info_classif z pakietu sklearn.feature_selection.feature_importances_.RFE). W tym algorytmie trenujemy klasyfikator na wszystkich cechach, wyrzucamy najsłabszą, trenujemy znowu i tak dalej.Typowo metody filter są najszybsze, ale dają najsłabszy wynik, natomiast metody wrapper są najwolniejsze i dają najlepszy wynik. Metody embedded są gdzieś pośrodku.
Dla zainteresowanych, inne znane i bardzo dobre algorytmy:
ReBATE): Wikipedia), artykuł "Benchmarking Relief-Based Feature Selection Methods"boruta_py): link 1, link 2from sklearn.feature_selection import RFE, SelectPercentile, mutual_info_classif
print(f"Without data selection\n")
print(f"Random forest:")
random_forest = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest.fit(X_reb, y_reb)
get_auroc_score(random_forest, X_test, y_test)
get_classification_report(
random_forest, X_test, y_test, "Classification report for Random forest:"
)
print(f"LGBM:")
LGBM = LGBMClassifier(importance_type="gain", random_state=0, n_jobs=-1)
LGBM.fit(X_reb, y_reb)
get_auroc_score(LGBM, X_test, y_test)
get_classification_report(LGBM, X_test, y_test, "Classification report for LGBM:")
Without data selection
Random forest:
AUROC score: 0.9047644274917003
Classification report for Random forest:
precision recall f1-score support
0 0.97 0.98 0.98 2002
1 0.53 0.40 0.46 99
accuracy 0.96 2101
macro avg 0.75 0.69 0.72 2101
weighted avg 0.95 0.96 0.95 2101
LGBM:
AUROC score: 0.9433748070111706
Classification report for LGBM:
precision recall f1-score support
0 0.98 0.98 0.98 2002
1 0.60 0.60 0.60 99
accuracy 0.96 2101
macro avg 0.79 0.79 0.79 2101
weighted avg 0.96 0.96 0.96 2101
print(f"Filter Methods\n")
filter_method = SelectPercentile(score_func=mutual_info_classif, percentile=80)
X_train_selected = filter_method.fit_transform(X_reb, y_reb)
X_test_selected = filter_method.transform(X_test)
print(f"Random forest:")
random_forest = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest.fit(X_train_selected, y_reb)
get_auroc_score(random_forest, X_test_selected, y_test)
get_classification_report(
random_forest, X_test_selected, y_test, "Classification report for Random forest:"
)
print(f"LGBM:")
LGBM = LGBMClassifier(importance_type="gain", random_state=0, n_jobs=-1)
LGBM.fit(X_train_selected, y_reb)
get_auroc_score(LGBM, X_test_selected, y_test)
get_classification_report(
LGBM, X_test_selected, y_test, "Classification report for LGBM:"
)
Filter Methods
Random forest:
AUROC score: 0.8951906679179406
Classification report for Random forest:
precision recall f1-score support
0 0.97 0.98 0.97 2002
1 0.46 0.42 0.44 99
accuracy 0.95 2101
macro avg 0.71 0.70 0.71 2101
weighted avg 0.95 0.95 0.95 2101
LGBM:
AUROC score: 0.9245602881966518
Classification report for LGBM:
precision recall f1-score support
0 0.98 0.97 0.98 2002
1 0.50 0.55 0.52 99
accuracy 0.95 2101
macro avg 0.74 0.76 0.75 2101
weighted avg 0.95 0.95 0.95 2101
# function for selecting 80% of the most importatnt features
def get_selected_features(feature_importances, X_train, X_test):
# get top 80% feature count
count = int(feature_importances.size * 0.8)
# calculate the threshold and cut the least important 20%
threshold = np.min(np.sort(feature_importances)[count:])
selected = feature_importances >= threshold
X_train_selected = X_train[:, selected]
X_test_selected = X_test[:, selected]
return X_train_selected, X_test_selected
print(f"Embedded Methods\n")
print(f"Random forest:")
X_train_forest, X_test_forest = get_selected_features(
forest_feature_importances, X_reb, X_test
)
random_forest = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest.fit(X_train_forest, y_reb)
get_auroc_score(random_forest, X_test_forest, y_test)
get_classification_report(
random_forest, X_test_forest, y_test, "Classification report for Random forest:"
)
print(f"LGBM:")
X_train_LGBM, X_test_LGBM = get_selected_features(
LGBM_feature_importances, X_reb, X_test
)
LGBM = LGBMClassifier(importance_type="gain", random_state=0, n_jobs=-1)
LGBM.fit(X_train_LGBM, y_reb)
get_auroc_score(LGBM, X_test_LGBM, y_test)
get_classification_report(LGBM, X_test_LGBM, y_test, "Classification report for LGBM:")
Embedded Methods
Random forest:
AUROC score: 0.9151757333575514
Classification report for Random forest:
precision recall f1-score support
0 0.98 0.97 0.97 2002
1 0.44 0.54 0.48 99
accuracy 0.95 2101
macro avg 0.71 0.75 0.73 2101
weighted avg 0.95 0.95 0.95 2101
LGBM:
AUROC score: 0.9339397975761612
Classification report for LGBM:
precision recall f1-score support
0 0.98 0.96 0.97 2002
1 0.45 0.62 0.52 99
accuracy 0.95 2101
macro avg 0.71 0.79 0.75 2101
weighted avg 0.96 0.95 0.95 2101
print(f"Wrapper Methods\n")
print(f"Random forest:")
random_forest = RandomForestClassifier(
n_estimators=500, criterion="entropy", random_state=0, n_jobs=-1
)
random_forest_RFE = RFE(estimator=random_forest, n_features_to_select=0.8)
random_forest_RFE.fit(X_reb, y_reb)
get_auroc_score(random_forest_RFE, X_test, y_test)
get_classification_report(
random_forest_RFE, X_test, y_test, "Classification report for Random forest:"
)
print(f"LGBM:")
LGBM = LGBMClassifier(importance_type="gain", random_state=0, n_jobs=-1)
LGBM_RFE = RFE(estimator=LGBM, n_features_to_select=0.8)
LGBM_RFE.fit(X_reb, y_reb)
get_auroc_score(LGBM_RFE, X_test, y_test)
get_classification_report(LGBM_RFE, X_test, y_test, "Classification report for LGBM:")
Wrapper Methods
Random forest:
AUROC score: 0.9105187741551377
Classification report for Random forest:
precision recall f1-score support
0 0.97 0.98 0.98 2002
1 0.51 0.41 0.46 99
accuracy 0.95 2101
macro avg 0.74 0.70 0.72 2101
weighted avg 0.95 0.95 0.95 2101
LGBM:
AUROC score: 0.9436876255058073
Classification report for LGBM:
precision recall f1-score support
0 0.98 0.98 0.98 2002
1 0.61 0.58 0.59 99
accuracy 0.96 2101
macro avg 0.79 0.78 0.79 2101
weighted avg 0.96 0.96 0.96 2101
Zgodnie z założeniami metody filter okazały się najsłabsze i tylko pogorszyły model. Metody Wrapper najbardziej usprawniły model, a metody embedded były gdzieś pośrodku. W przypadku innych metryk: precyzja i pokrycie, selekcja cech na ogół negatywnie na nie wpłynęła, bądź lekko polepszyła.